home *** CD-ROM | disk | FTP | other *** search
/ Visual Basic Graphics Programming (2nd Edition) / Visual Basic Graphics Programming 2nd Edition.iso / Src / Ch17 / Matrix.bas < prev    next >
BASIC Source File  |  1999-06-05  |  17KB  |  580 lines

  1. Attribute VB_Name = "Matrix"
  2. Option Explicit
  3.  
  4. ' Routines for making and manipulating matrices.
  5.  
  6. Public Const PI = 3.14159265358979
  7. Public Const INFINITY = 2147483647
  8.  
  9. Public Enum ProjectionTypes
  10.     project_Parallel
  11.     project_Perspective
  12. End Enum
  13.  
  14. Public Type Point3D
  15.     Coord(1 To 4) As Single
  16.     Trans(1 To 4) As Single
  17. End Type
  18.  
  19. Public Type Segment3D
  20.     pt1 As Integer
  21.     pt2 As Integer
  22. End Type
  23. ' Convert the spherical coordinates into
  24. ' Cartesian coordinates.
  25. Public Sub m3SphericalToCartesian(ByVal r As Single, ByVal theta As Single, ByVal phi As Single, X As Single, Y As Single, Z As Single)
  26. Dim r2 As Single
  27.     
  28.     ' Create a line to the center of projection.
  29.     Y = r * Sin(phi)
  30.     r2 = r * Cos(phi)
  31.     X = r2 * Cos(theta)
  32.     Z = r2 * Sin(theta)
  33. End Sub
  34. ' Create a transformation matrix for an oblique
  35. ' projection onto the X-Y plane.
  36. Public Sub m3ObliqueXY(M() As Single, ByVal S As Single, ByVal theta As Single)
  37.     m3Identity M
  38.     M(3, 1) = -S * Cos(theta)
  39.     M(3, 2) = -S * Sin(theta)
  40.     M(3, 3) = 0
  41. End Sub
  42.  
  43.  
  44. ' Create a transformation matrix for orthographic
  45. ' projection along the X axis.
  46. Public Sub m3OrthoSide(M() As Single)
  47.     m3Identity M
  48.     M(1, 1) = 0
  49.     M(3, 1) = -1
  50.     M(3, 3) = 0
  51. End Sub
  52. ' Create a transformation matrix for orthographic
  53. ' projection along the Y axis.
  54. Public Sub m3OrthoTop(M() As Single)
  55.     m3Identity M
  56.     M(2, 2) = 0
  57.     M(3, 2) = -1
  58.     M(3, 3) = 0
  59. End Sub
  60.  
  61. ' Create a transformation matrix for orthographic
  62. ' projection along the Z axis.
  63. Public Sub m3OrthoFront(M() As Single)
  64.     m3Identity M
  65.     M(3, 3) = 0
  66. End Sub
  67.  
  68. ' Create an identity matrix.
  69. Public Sub m3Identity(M() As Single)
  70. Dim i As Integer
  71. Dim j As Integer
  72.  
  73.     For i = 1 To 4
  74.         For j = 1 To 4
  75.             If i = j Then
  76.                 M(i, j) = 1
  77.             Else
  78.                 M(i, j) = 0
  79.             End If
  80.         Next j
  81.     Next i
  82. End Sub
  83.  
  84. ' Normalize a 3-D point vector.
  85. Public Sub m3NormalizeCoords(X As Single, Y As Single, Z As Single, S As Single)
  86.     X = X / S
  87.     Y = Y / S
  88.     Z = Z / S
  89.     S = 1
  90. End Sub
  91.  
  92. ' Normalize a 3-D point vector.
  93. Public Sub m3NormalizePoint(P() As Single)
  94. Dim i As Integer
  95. Dim value As Single
  96.  
  97.     value = P(4)
  98.     For i = 1 To 3
  99.         P(i) = P(i) / value
  100.     Next i
  101.     P(4) = 1
  102. End Sub
  103.  
  104.  
  105. ' Normalize a 3-D transformation matrix.
  106. Public Sub m3NormalizeMatrix(M() As Single)
  107. Dim i As Integer
  108. Dim j As Integer
  109. Dim value As Single
  110.  
  111.     value = M(4, 4)
  112.     For i = 1 To 4
  113.         For j = 1 To 4
  114.             M(i, j) = M(i, j) / value
  115.         Next j
  116.     Next i
  117. End Sub
  118.  
  119.  
  120.  
  121.  
  122. ' Create a 3-D transformation matrix for a
  123. ' perspective projection along the Z axis onto
  124. ' the X-Y plane with focus at the origin and the
  125. ' center of projection at distance (0, 0, D).
  126. Public Sub project_PerspectiveXY(M() As Single, ByVal D As Single)
  127.     m3Identity M
  128.     If D <> 0 Then M(3, 4) = -1 / D
  129. End Sub
  130.  
  131. ' Create a 3-D transformation matrix for a
  132. ' projection with:
  133. '       center of projection    (cx, cy, cz)
  134. '       focus                   (fx, fy, fx)
  135. '       UP vector               <ux, yx, uz>
  136. ' ptype should be project_Perspective or project_Parallel.
  137. Public Sub m3Project(M() As Single, ByVal ptype As ProjectionTypes, ByVal Cx As Single, ByVal Cy As Single, ByVal Cz As Single, ByVal Fx As Single, ByVal Fy As Single, ByVal Fz As Single, ByVal ux As Single, ByVal uy As Single, ByVal uz As Single)
  138. Static M1(1 To 4, 1 To 4) As Single
  139. Static M2(1 To 4, 1 To 4) As Single
  140. Static M3(1 To 4, 1 To 4) As Single
  141. Static M4(1 To 4, 1 To 4) As Single
  142. Static M5(1 To 4, 1 To 4) As Single
  143. Static M12(1 To 4, 1 To 4) As Single
  144. Static M34(1 To 4, 1 To 4) As Single
  145. Static M1234(1 To 4, 1 To 4) As Single
  146. Dim sin1 As Single
  147. Dim cos1 As Single
  148. Dim sin2 As Single
  149. Dim cos2 As Single
  150. Dim sin3 As Single
  151. Dim cos3 As Single
  152. Dim A As Single
  153. Dim B As Single
  154. Dim C As Single
  155. Dim d1 As Single
  156. Dim d2 As Single
  157. Dim d3 As Single
  158. Dim up1(1 To 4) As Single
  159. Dim up2(1 To 4) As Single
  160.  
  161.     ' Translate the focus to the origin.
  162.     m3Translate M1, -Fx, -Fy, -Fz
  163.  
  164.     A = Cx - Fx
  165.     B = Cy - Fy
  166.     C = Cz - Fz
  167.     d1 = Sqr(A * A + C * C)
  168.     If d1 <> 0 Then
  169.         sin1 = -A / d1
  170.         cos1 = C / d1
  171.     End If
  172.     d2 = Sqr(A * A + B * B + C * C)
  173.     If d2 <> 0 Then
  174.         sin2 = B / d2
  175.         cos2 = d1 / d2
  176.     End If
  177.     
  178.     ' Rotate around the Y axis to place the
  179.     ' center of projection in the Y-Z plane.
  180.     m3Identity M2
  181.     
  182.     ' If d1 = 0 then the center of projection
  183.     ' already lies in the Y axis and thus the Y-Z plane.
  184.     If d1 <> 0 Then
  185.         M2(1, 1) = cos1
  186.         M2(1, 3) = -sin1
  187.         M2(3, 1) = sin1
  188.         M2(3, 3) = cos1
  189.     End If
  190.     
  191.     ' Rotate around the X axis to place the
  192.     ' center of projection in the Z axis.
  193.     m3Identity M3
  194.     
  195.     ' If d2 = 0 then the center of projection
  196.     ' lies at the origin. This makes projection
  197.     ' impossible.
  198.     If d2 <> 0 Then
  199.         M3(2, 2) = cos2
  200.         M3(2, 3) = sin2
  201.         M3(3, 2) = -sin2
  202.         M3(3, 3) = cos2
  203.     End If
  204.     
  205.     ' Apply the rotations to the UP vector.
  206.     up1(1) = ux
  207.     up1(2) = uy
  208.     up1(3) = uz
  209.     up1(4) = 1
  210.     m3Apply up1, M2, up2
  211.     m3Apply up2, M3, up1
  212.  
  213.     ' Rotate around the Z axis to put the UP
  214.     ' vector in the Y-Z plane.
  215.     d3 = Sqr(up1(1) * up1(1) + up1(2) * up1(2))
  216.     m3Identity M4
  217.     
  218.     ' If d3 = 0 then the UP vector is a zero
  219.     ' vector so do nothing.
  220.     If d3 <> 0 Then
  221.         sin3 = up1(1) / d3
  222.         cos3 = up1(2) / d3
  223.         M4(1, 1) = cos3
  224.         M4(1, 2) = sin3
  225.         M4(2, 1) = -sin3
  226.         M4(2, 2) = cos3
  227.     End If
  228.     
  229.     ' Project.
  230.     If ptype = project_Perspective And d2 <> 0 Then
  231.         project_PerspectiveXY M5, d2
  232.     Else
  233.         m3Identity M5
  234.     End If
  235.  
  236.     ' Combine the transformations.
  237.     m3MatMultiply M12, M1, M2
  238.     m3MatMultiply M34, M3, M4
  239.     m3MatMultiply M1234, M12, M34
  240.     If ptype = project_Perspective Then
  241.         m3MatMultiplyFull M, M1234, M5
  242.     Else
  243.         m3MatMultiply M, M1234, M5
  244.     End If
  245. End Sub
  246.  
  247.  
  248.  
  249. ' Create a 3-D transformation matrix for a
  250. ' perspective projection with:
  251. '       center of projection    (r, phi, theta)
  252. '       focus                   (fx, fy, fx)
  253. '       up vector               <ux, uy, uz>
  254. ' ptype should be project_Perspective or project_Parallel.
  255. Public Sub m3PProject(M() As Single, ByVal ptype As ProjectionTypes, ByVal r As Single, ByVal phi As Single, ByVal theta As Single, ByVal Fx As Single, ByVal Fy As Single, ByVal Fz As Single, ByVal ux As Single, ByVal uy As Single, ByVal uz As Single)
  256. Dim Cx As Single
  257. Dim Cy As Single
  258. Dim Cz As Single
  259. Dim r2 As Single
  260.  
  261.     ' Convert to Cartesian coordinates.
  262.     Cy = r * Sin(phi)
  263.     r2 = r * Cos(phi)
  264.     Cx = r2 * Cos(theta)
  265.     Cz = r2 * Sin(theta)
  266.     m3Project M, ptype, Cx, Cy, Cz, Fx, Fy, Fz, ux, uy, uz
  267. End Sub
  268.  
  269. ' Create a transformation matrix for reflecting
  270. ' across the plane passing through (p1, p2, p3)
  271. ' with normal vector <n1, n2, n3>.
  272. Public Sub m3Reflect(M() As Single, ByVal p1 As Single, ByVal p2 As Single, ByVal p3 As Single, ByVal n1 As Single, ByVal n2 As Single, ByVal n3 As Single)
  273. Dim T(1 To 4, 1 To 4) As Single     ' Translate.
  274. Dim R1(1 To 4, 1 To 4) As Single    ' Rotate 1.
  275. Dim r2(1 To 4, 1 To 4) As Single    ' Rotate 2.
  276. Dim S(1 To 4, 1 To 4) As Single     ' Reflect.
  277. Dim R2i(1 To 4, 1 To 4) As Single   ' Unrotate 2.
  278. Dim R1i(1 To 4, 1 To 4) As Single   ' Unrotate 1.
  279. Dim Ti(1 To 4, 1 To 4) As Single    ' Untranslate.
  280. Dim D As Single
  281. Dim l As Single
  282. Dim M12(1 To 4, 1 To 4) As Single
  283. Dim M34(1 To 4, 1 To 4) As Single
  284. Dim M1234(1 To 4, 1 To 4) As Single
  285. Dim M56(1 To 4, 1 To 4) As Single
  286. Dim M567(1 To 4, 1 To 4) As Single
  287.  
  288.     ' Translate the plane to the origin.
  289.     m3Translate T, -p1, -p2, -p3
  290.     m3Translate Ti, p1, p2, p3
  291.  
  292.     ' Rotate around Z-axis until the normal is in
  293.     ' the Y-Z plane.
  294.     m3Identity R1
  295.     D = Sqr(n1 * n1 + n2 * n2)
  296.     R1(1, 1) = n2 / D
  297.     R1(1, 2) = n1 / D
  298.     R1(2, 1) = -R1(1, 2)
  299.     R1(2, 2) = R1(1, 1)
  300.     
  301.     m3Identity R1i
  302.     R1i(1, 1) = R1(1, 1)
  303.     R1i(1, 2) = -R1(1, 2)
  304.     R1i(2, 1) = -R1(2, 1)
  305.     R1i(2, 2) = R1(2, 2)
  306.     
  307.     ' Rotate around the X-axis until the normal
  308.     ' lies along the Y axis.
  309.     m3Identity r2
  310.     l = Sqr(n1 * n1 + n2 * n2 + n3 * n3)
  311.     r2(2, 2) = D / l
  312.     r2(2, 3) = -n3 / l
  313.     r2(3, 2) = -r2(2, 3)
  314.     r2(3, 3) = r2(2, 2)
  315.     
  316.     m3Identity R2i
  317.     R2i(2, 2) = r2(2, 2)
  318.     R2i(2, 3) = -r2(2, 3)
  319.     R2i(3, 2) = -r2(3, 2)
  320.     R2i(3, 3) = r2(3, 3)
  321.  
  322.     ' Reflect across the X-Z plane.
  323.     m3Identity S
  324.     S(2, 2) = -1
  325.  
  326.     ' Combine the matrices.
  327.     m3MatMultiply M12, T, R1
  328.     m3MatMultiply M34, r2, S
  329.     m3MatMultiply M1234, M12, M34
  330.     m3MatMultiply M56, R2i, R1i
  331.     m3MatMultiply M567, M56, Ti
  332.     m3MatMultiply M, M1234, M567
  333. End Sub
  334.  
  335.  
  336. ' Create a transformation atrix for rotating
  337. ' through angle theta around a line passing
  338. ' through (p1, p2, p3) in direction <d1, d2, d3>.
  339. ' Theta is measured counterclockwise as you look
  340. ' down the line opposite the line's direction.
  341. Public Sub m3LineRotate(M() As Single, ByVal p1 As Single, ByVal p2 As Single, ByVal p3 As Single, ByVal d1 As Single, ByVal d2 As Single, ByVal d3 As Single, ByVal theta As Single)
  342. Dim T(1 To 4, 1 To 4) As Single     ' Translate.
  343. Dim R1(1 To 4, 1 To 4) As Single    ' Rotate 1.
  344. Dim r2(1 To 4, 1 To 4) As Single    ' Rotate 2.
  345. Dim Rot3(1 To 4, 1 To 4) As Single  ' Rotate.
  346. Dim R2i(1 To 4, 1 To 4) As Single   ' Unrotate 2.
  347. Dim R1i(1 To 4, 1 To 4) As Single   ' Unrotate 1.
  348. Dim Ti(1 To 4, 1 To 4) As Single    ' Untranslate.
  349. Dim D As Single
  350. Dim l As Single
  351. Dim M12(1 To 4, 1 To 4) As Single
  352. Dim M34(1 To 4, 1 To 4) As Single
  353. Dim M1234(1 To 4, 1 To 4) As Single
  354. Dim M56(1 To 4, 1 To 4) As Single
  355. Dim M567(1 To 4, 1 To 4) As Single
  356.  
  357.     ' Translate the plane to the origin.
  358.     m3Translate T, -p1, -p2, -p3
  359.     m3Translate Ti, p1, p2, p3
  360.  
  361.     ' Rotate around Z-axis until the line is in
  362.     ' the Y-Z plane.
  363.     m3Identity R1
  364.     D = Sqr(d1 * d1 + d2 * d2)
  365.     R1(1, 1) = d2 / D
  366.     R1(1, 2) = d1 / D
  367.     R1(2, 1) = -R1(1, 2)
  368.     R1(2, 2) = R1(1, 1)
  369.     
  370.     m3Identity R1i
  371.     R1i(1, 1) = R1(1, 1)
  372.     R1i(1, 2) = -R1(1, 2)
  373.     R1i(2, 1) = -R1(2, 1)
  374.     R1i(2, 2) = R1(2, 2)
  375.     
  376.     ' Rotate around the X-axis until the line
  377.     ' lies along the Y axis.
  378.     m3Identity r2
  379.     l = Sqr(d1 * d1 + d2 * d2 + d3 * d3)
  380.     r2(2, 2) = D / l
  381.     r2(2, 3) = -d3 / l
  382.     r2(3, 2) = -r2(2, 3)
  383.     r2(3, 3) = r2(2, 2)
  384.     
  385.     m3Identity R2i
  386.     R2i(2, 2) = r2(2, 2)
  387.     R2i(2, 3) = -r2(2, 3)
  388.     R2i(3, 2) = -r2(3, 2)
  389.     R2i(3, 3) = r2(3, 3)
  390.  
  391.     ' Rotate around the line (Y axis).
  392.     m3YRotate Rot3, theta
  393.  
  394.     ' Combine the matrices.
  395.     m3MatMultiply M12, T, R1
  396.     m3MatMultiply M34, r2, Rot3
  397.     m3MatMultiply M1234, M12, M34
  398.     m3MatMultiply M56, R2i, R1i
  399.     m3MatMultiply M567, M56, Ti
  400.     m3MatMultiply M, M1234, M567
  401. End Sub
  402.  
  403. ' Create a 3-D transformation matrix for scaling
  404. ' by scale factors Sx, Sy, and Sz.
  405. Public Sub m3Scale(M() As Single, ByVal Sx As Single, ByVal Sy As Single, ByVal Sz As Single)
  406.     m3Identity M
  407.     M(1, 1) = Sx
  408.     M(2, 2) = Sy
  409.     M(3, 3) = Sz
  410. End Sub
  411.  
  412. ' Create a 3-D transformation matrix for
  413. ' translation by Tx, Ty, and Tz.
  414. Public Sub m3Translate(M() As Single, ByVal Tx As Single, ByVal Ty As Single, ByVal Tz As Single)
  415.     m3Identity M
  416.     M(4, 1) = Tx
  417.     M(4, 2) = Ty
  418.     M(4, 3) = Tz
  419. End Sub
  420.  
  421. ' Create a 3-D transformation matrix for rotation
  422. ' around the X axis (angle measured in radians).
  423. Public Sub m3XRotate(M() As Single, ByVal theta As Single)
  424.     m3Identity M
  425.     M(2, 2) = Cos(theta)
  426.     M(3, 3) = M(2, 2)
  427.     M(2, 3) = Sin(theta)
  428.     M(3, 2) = -M(2, 3)
  429. End Sub
  430.  
  431. ' Create a 3-D transformation matrix for rotation
  432. ' around the Y axis (angle measured in radians).
  433. Public Sub m3YRotate(M() As Single, ByVal theta As Single)
  434.     m3Identity M
  435.     M(1, 1) = Cos(theta)
  436.     M(3, 3) = M(1, 1)
  437.     M(3, 1) = Sin(theta)
  438.     M(1, 3) = -M(3, 1)
  439. End Sub
  440.  
  441. ' Create a 3-D transformation matrix for rotation
  442. ' around the Z axis (angle measured in radians).
  443. Public Sub m3ZRotate(M() As Single, ByVal theta As Single)
  444.     m3Identity M
  445.     M(1, 1) = Cos(theta)
  446.     M(2, 2) = M(1, 1)
  447.     M(1, 2) = Sin(theta)
  448.     M(2, 1) = -M(1, 2)
  449. End Sub
  450.  
  451. ' Create a matrix that rotates around the Y axis
  452. ' so the point (x, y, z) lies in the X-Z plane.
  453. Public Sub m3YRotateIntoXZ(Result() As Single, ByVal X As Single, ByVal Y As Single, ByVal Z As Single)
  454. Dim D As Single
  455.  
  456.     m3Identity Result
  457.     D = Sqr(X * X + Y * Y)
  458.     Result(1, 1) = X / D
  459.     Result(1, 2) = -Y / D
  460.     Result(2, 1) = -Result(1, 2)
  461.     Result(2, 2) = Result(1, 1)
  462. End Sub
  463.  
  464. ' Set copy = orig.
  465. Public Sub m3MatCopy(copy() As Single, orig() As Single)
  466. Dim i As Integer
  467. Dim j As Integer
  468.  
  469.     For i = 1 To 4
  470.         For j = 1 To 4
  471.             copy(i, j) = orig(i, j)
  472.         Next j
  473.     Next i
  474. End Sub
  475.  
  476. ' Apply a transformation matrix to a point where
  477. ' the transformation may not have 0, 0, 0, 1 in
  478. ' its final column. Normalize only the X and Y
  479. ' components of the result to preserve the Z
  480. ' information.
  481. Public Sub m3ApplyFull(V() As Single, M() As Single, Result() As Single)
  482. Dim i As Integer
  483. Dim j As Integer
  484. Dim value As Single
  485.  
  486.     For i = 1 To 4
  487.         value = 0#
  488.         For j = 1 To 4
  489.             value = value + V(j) * M(j, i)
  490.         Next j
  491.         Result(i) = value
  492.     Next i
  493.  
  494.     ' Renormalize the point.
  495.     ' Note that value still holds Result(4).
  496.     If value <> 0 Then
  497.         Result(1) = Result(1) / value
  498.         Result(2) = Result(2) / value
  499.         ' Do not transform the Z component.
  500.     Else
  501.         ' Make the Z value greater than that of
  502.         ' the center of projection so the point
  503.         ' will be clipped.
  504.         Result(3) = INFINITY
  505.     End If
  506.     Result(4) = 1#
  507. End Sub
  508.  
  509.  
  510.  
  511.  
  512. ' Apply a transformation matrix to a point.
  513. Public Sub m3Apply(V() As Single, M() As Single, Result() As Single)
  514.     Result(1) = V(1) * M(1, 1) + _
  515.                 V(2) * M(2, 1) + _
  516.                 V(3) * M(3, 1) + M(4, 1)
  517.     Result(2) = V(1) * M(1, 2) + _
  518.                 V(2) * M(2, 2) + _
  519.                 V(3) * M(3, 2) + M(4, 2)
  520.     Result(3) = V(1) * M(1, 3) + _
  521.                 V(2) * M(2, 3) + _
  522.                 V(3) * M(3, 3) + M(4, 3)
  523.     Result(4) = 1#
  524. End Sub
  525.  
  526. ' Multiply two matrices together. The matrices
  527. ' may not contain 0, 0, 0, 1 in their last
  528. ' columns.
  529. Public Sub m3MatMultiplyFull(Result() As Single, A() As Single, B() As Single)
  530. Dim i As Integer
  531. Dim j As Integer
  532. Dim k As Integer
  533. Dim value As Single
  534.  
  535.     For i = 1 To 4
  536.         For j = 1 To 4
  537.             value = 0#
  538.             For k = 1 To 4
  539.                 value = value + A(i, k) * B(k, j)
  540.             Next k
  541.             Result(i, j) = value
  542.         Next j
  543.     Next i
  544. End Sub
  545. ' Multiply two matrices together.
  546. Public Sub m3MatMultiply(Result() As Single, A() As Single, B() As Single)
  547.     Result(1, 1) = A(1, 1) * B(1, 1) + A(1, 2) * B(2, 1) + A(1, 3) * B(3, 1)
  548.     Result(1, 2) = A(1, 1) * B(1, 2) + A(1, 2) * B(2, 2) + A(1, 3) * B(3, 2)
  549.     Result(1, 3) = A(1, 1) * B(1, 3) + A(1, 2) * B(2, 3) + A(1, 3) * B(3, 3)
  550.     Result(1, 4) = 0#
  551.     Result(2, 1) = A(2, 1) * B(1, 1) + A(2, 2) * B(2, 1) + A(2, 3) * B(3, 1)
  552.     Result(2, 2) = A(2, 1) * B(1, 2) + A(2, 2) * B(2, 2) + A(2, 3) * B(3, 2)
  553.     Result(2, 3) = A(2, 1) * B(1, 3) + A(2, 2) * B(2, 3) + A(2, 3) * B(3, 3)
  554.     Result(2, 4) = 0#
  555.     Result(3, 1) = A(3, 1) * B(1, 1) + A(3, 2) * B(2, 1) + A(3, 3) * B(3, 1)
  556.     Result(3, 2) = A(3, 1) * B(1, 2) + A(3, 2) * B(2, 2) + A(3, 3) * B(3, 2)
  557.     Result(3, 3) = A(3, 1) * B(1, 3) + A(3, 2) * B(2, 3) + A(3, 3) * B(3, 3)
  558.     Result(3, 4) = 0#
  559.     Result(4, 1) = A(4, 1) * B(1, 1) + A(4, 2) * B(2, 1) + A(4, 3) * B(3, 1) + B(4, 1)
  560.     Result(4, 2) = A(4, 1) * B(1, 2) + A(4, 2) * B(2, 2) + A(4, 3) * B(3, 2) + B(4, 2)
  561.     Result(4, 3) = A(4, 1) * B(1, 3) + A(4, 2) * B(2, 3) + A(4, 3) * B(3, 3) + B(4, 3)
  562.     Result(4, 4) = 1#
  563. End Sub
  564.  
  565. ' Compute the cross product of two vectors.
  566. ' Set <x, y, z> = <x1, y1, z1> X <x2, y2, z2>.
  567. Public Sub m3Cross(X As Single, Y As Single, Z As Single, ByVal x1 As Single, ByVal y1 As Single, ByVal z1 As Single, ByVal x2 As Single, ByVal y2 As Single, ByVal z2 As Single)
  568.     X = y1 * z2 - z1 * y2
  569.     Y = z1 * x2 - x1 * z2
  570.     Z = x1 * y2 - y1 * x2
  571. End Sub
  572.  
  573. ' Give the vector the indicated length.
  574. Public Sub m3SizeVector(ByVal l As Single, X As Single, Y As Single, Z As Single)
  575.     l = l / Sqr(X * X + Y * Y + Z * Z)
  576.     X = X * l
  577.     Y = Y * l
  578.     Z = Z * l
  579. End Sub
  580.